#ifndef _MESH_
#define _MESH_

#include "Domain.h"
#include "Dofs.h"
#include "Grid.h"
#include <map>
#include <functional>
#include <fstream>
using namespace std;
/**
 * @brief 定义域划分成网格Mesh,Mesh包含小网格Grid,每个Grid的类型由TemplateElement确定
 *
 * @tparam int DIM 维度
 */
template <unsigned int DIM>
class Mesh
{
public:
    /**
     * @brief 默认构造函数
     */
    Mesh() = default;

    /**
     * @brief 默认析构函数
     */
    virtual ~Mesh() = default;
    /**
     * @brief 返回自由度的总个数
     * @details _GridList.size()*getNumDofs();
     * @return
     */
    virtual int getTotalNumDofs() = 0;
    /**
     * @brief 返回网格的总个数
     * @details _GridList.size()
     * @return [description]
     */
    virtual int getNumGrid() = 0;
    /**
     * @brief 返回存储Grid的vector
     * @details _GridList
     * @return [description]
     */
    virtual vector<Grid<DIM>> &getGrid() = 0;
    /**
     * @brief 返回序号为_index的Grid
     *
     * @param _index Grid的序号
     * @return _GridList[i];
     */
    virtual Grid<DIM> &getGrid(int _index) = 0;
    /**
     * @brief   返回vector，存储边界上的自由度
     * @details 对所有grid,依据seg找出边界的grid,再确定对应边界的自由度下标是哪几个
     *
     * @param component 试图对边界分段，并标号
     * @return 存储边界上自由度的vector
     */
    virtual vector<Grid<2>> getBoundaryGrid() = 0;
    virtual void setBoundaryIndex() = 0;
    virtual void distributeDofs() = 0;
    virtual vector<pair<int, int>> getIndexofGrid(int i) = 0;
    virtual vector<int> &getSegment() = 0;
};

class Q1Mesh : public Mesh<2>
{
public:
    Q1Mesh() = default;
    /**
     * @brief 含参构造函数,根据给定的定义域和划分段数来划分网格。
     * @param domain 定义域
     * @param seg 每个维度要划分的段数
     */
    Q1Mesh(Domain<2> *_domain, vector<int> _seg);
    ~Q1Mesh() = default;

    int getTotalNumDofs();
    int getNumGrid();
    vector<Grid<2>> &getGrid();
    Grid<2> &getGrid(int i);
    vector<Grid<2>> getBoundaryGrid();
    vector<int> &getSegment();
    void distributeDofs();
    void setBoundaryIndex();
    vector<pair<int, int>> getIndexofGrid(int i);

protected:
    vector<Grid<2>> gridList;           /*< 存放所有的Grid*/
    vector<pair<int, int>> IndexofGrid; /*< 存放grid上的所有index */
    vector<int> seg;                    /*< 存放各维度分段数*/
};

Q1Mesh::Q1Mesh(Domain<2> *_domain, vector<int> _seg)
{
    assert(_seg.size() == 2);
    Point<2> p0 = _domain->getVertex(0);
    Point<2> p2 = _domain->getVertex(2);
    double xLen = (p2[0] - p0[0]) / (1.0 * _seg[0]);
    double yLen = (p2[1] - p0[1]) / (1.0 * _seg[1]);
    for (int j = 0; j < _seg[1]; ++j)
    {
        for (int i = 0; i < _seg[0]; ++i)
        {
            vector<Point<2>> pnts;
            pnts.push_back({p0[0] + i * xLen, p0[1] + j * yLen});
            pnts.push_back({p0[0] + (i + 1) * xLen, p0[1] + j * yLen});
            pnts.push_back({p0[0] + (i + 1) * xLen, p0[1] + (j + 1) * yLen});
            pnts.push_back({p0[0] + i * xLen, p0[1] + (j + 1) * yLen});
            Grid<2> grid(pnts);
            gridList.push_back(grid);
        }
    }
    seg = _seg;
}

int Q1Mesh::getNumGrid()
{
    return gridList.size();
}

Grid<2> &Q1Mesh::getGrid(int i)
{
    return gridList[i];
}

vector<Grid<2>> &Q1Mesh::getGrid()
{
    return gridList;
}

vector<Grid<2>> Q1Mesh::getBoundaryGrid()
{
    vector<Grid<2>> boundary_grid;
    int k = 0;
    for (int i = 0; i < seg[0] - 1; i++)
    {
        boundary_grid.push_back(gridList[k]);
        k++;
    }
    k = seg[0];
    for (int i = 0; i < seg[1] - 1; i++)
    {
        boundary_grid.push_back(gridList[k]);
        k += seg[0];
    }
    k = seg[0] - 1;
    for (int i = 0; i < seg[1] - 1; i++)
    {
        boundary_grid.push_back(gridList[k]);
        k += seg[0];
    }
    k = seg[0] * (seg[1] - 1);
    for (int i = 0; i < seg[0] - 1; i++)
    {
        boundary_grid.push_back(gridList[k]);
        k++;
    }
    return boundary_grid;
}

vector<int> &Q1Mesh::getSegment()
{
    return seg;
}

vector<pair<int, int>> Q1Mesh::getIndexofGrid(int i)
{
    vector<pair<int, int>> idx;
    int a = seg[0];
    idx.push_back(IndexofGrid[i + i / a]);
    idx.push_back(IndexofGrid[i + i / a + 1]);
    idx.push_back(IndexofGrid[i + i / a + a + 2]);
    idx.push_back(IndexofGrid[i + i / a + a + 1]);
    return idx;
}

void Q1Mesh::setBoundaryIndex()
{
    for (int k = 0; k < seg[0]; k++)
        IndexofGrid[k].second = 1;
    for (int k = seg[1] + 1; k <= (seg[0] + 1) * seg[1]; k = k + seg[0] + 1)
        IndexofGrid[k].second = 1;
    for (int k = seg[0]; k < (seg[0] + 1) * (seg[1] + 1); k = k + seg[0] + 1)
        IndexofGrid[k].second = 1;
    for (int k = (seg[0] + 1) * seg[1] + 1; k < (seg[0] + 1) * (seg[1] + 1); k++)
        IndexofGrid[k].second = 1;
}

void Q1Mesh::distributeDofs()
{
    for (int i = 0; i < (seg[0] + 1) * (seg[1] + 1); i++)
    {
        IndexofGrid.push_back({i, 0});
    }
}

int Q1Mesh::getTotalNumDofs()
{
    return (seg[0] + 1) * (seg[1] + 1);
}

class Q2Mesh : public Mesh<2>
{
public:
    Q2Mesh() = default;
    /**
     * @brief 含参构造函数,根据给定的定义域和划分段数来划分网格。
     * @param domain 定义域
     * @param seg 每个维度要划分的段数
     */
    Q2Mesh(Domain<2> *_domain, vector<int> _seg);
    ~Q2Mesh() = default;

    int getTotalNumDofs();
    int getNumGrid();
    vector<Grid<2>> &getGrid();
    Grid<2> &getGrid(int i);
    vector<Grid<2>> getBoundaryGrid();
    vector<int> &getSegment();
    void distributeDofs();
    void setBoundaryIndex();
    vector<pair<int, int>> getIndexofGrid(int i);

protected:
    vector<Grid<2>> gridList;           /*< 存放所有的Grid*/
    vector<pair<int, int>> IndexofGrid; /*< 存放grid上的所有index */
    vector<int> seg;                    /*< 存放各维度分段数*/
};

Q2Mesh::Q2Mesh(Domain<2> *_domain, vector<int> _seg)
{
    assert(_seg.size() == 2);
    Point<2> p0 = _domain->getVertex(0);
    Point<2> p2 = _domain->getVertex(2);
    double xLen = (p2[0] - p0[0]) / (1.0 * _seg[0]);
    double yLen = (p2[1] - p0[1]) / (1.0 * _seg[1]);
    for (int j = 0; j < _seg[1]; ++j)
    {
        for (int i = 0; i < _seg[0]; ++i)
        {
            vector<Point<2>> pnts;
            pnts.push_back({p0[0] + i * xLen, p0[1] + j * yLen});
            pnts.push_back({p0[0] + (i + 1) * xLen, p0[1] + j * yLen});
            pnts.push_back({p0[0] + (i + 1) * xLen, p0[1] + (j + 1) * yLen});
            pnts.push_back({p0[0] + i * xLen, p0[1] + (j + 1) * yLen});
            Grid<2> grid(pnts);
            gridList.push_back(grid);
        }
    }
    seg = _seg;
}

int Q2Mesh::getNumGrid()
{
    return gridList.size();
}

Grid<2> &Q2Mesh::getGrid(int i)
{
    return gridList[i];
}

vector<Grid<2>> &Q2Mesh::getGrid()
{
    return gridList;
}

vector<int> &Q2Mesh::getSegment()
{
    return seg;
}

int Q2Mesh::getTotalNumDofs()
{
    return (2 * seg[0] + 1) * (2 * seg[1] + 1);
}

void Q2Mesh::distributeDofs()
{
    for (int i = 0; i < (2 * seg[0] + 1) * (2 * seg[1] + 1); i++)
    {
        IndexofGrid.push_back({i, 0});
    }
}

void Q2Mesh::setBoundaryIndex()
{
    int k = 0;
    for (int i = 0; i < 2 * seg[0]; i++)
    {
        IndexofGrid[k].second = 1;
        k++;
    }
    k = 2 * seg[0] + 1;
    for (int i = 0; i < 2 * seg[1]; i++)
    {
        IndexofGrid[k].second = 1;
        k += 2 * seg[0] + 1;
    }
    k = 2 * seg[0];
    for (int i = 0; i < 2 * seg[1]; i++)
    {
        IndexofGrid[k].second = 1;
        k += 2 * seg[0] + 1;
    }
    k = (2 * seg[0] + 1) * (2 * seg[1]) + 1;
    for (int i = 0; i < 2 * seg[0]; i++)
    {
        IndexofGrid[k].second = 1;
        k++;
    }
}

vector<pair<int, int>> Q2Mesh::getIndexofGrid(int i)
{
    vector<pair<int, int>> idx;
    int a = seg[0];
    idx.push_back(IndexofGrid[i / a * (4 * a + 2) + i % a * 2]);
    idx.push_back(IndexofGrid[i / a * (4 * a + 2) + i % a * 2 + 1]);
    idx.push_back(IndexofGrid[i / a * (4 * a + 2) + i % a * 2 + 2]);
    idx.push_back(IndexofGrid[i / a * (4 * a + 2) + 2 * a + 1 + i % a * 2 + 2]);
    idx.push_back(IndexofGrid[(i / a + 1) * (4 * a + 2) + i % a * 2 + 2]);
    idx.push_back(IndexofGrid[(i / a + 1) * (4 * a + 2) + i % a * 2 + 1]);
    idx.push_back(IndexofGrid[(i / a + 1) * (4 * a + 2) + i % a * 2]);
    idx.push_back(IndexofGrid[i / a * (4 * a + 2) + 2 * a + 1 + i % a * 2]);
    idx.push_back(IndexofGrid[i / a * (4 * a + 2) + 2 * a + 1 + i % a * 2 + 1]);
    return idx;
}

vector<Grid<2>> Q2Mesh::getBoundaryGrid()
{
    vector<Grid<2>> boundary_grid;
    int k = 0;
    for (int i = 0; i < seg[0] - 1; i++)
    {
        boundary_grid.push_back(gridList[k]);
        k++;
    }
    k = seg[0];
    for (int i = 0; i < seg[1] - 1; i++)
    {
        boundary_grid.push_back(gridList[k]);
        k += seg[0];
    }
    k = seg[0] - 1;
    for (int i = 0; i < seg[1] - 1; i++)
    {
        boundary_grid.push_back(gridList[k]);
        k += seg[0];
    }
    k = seg[0] * (seg[1] - 1);
    for (int i = 0; i < seg[0] - 1; i++)
    {
        boundary_grid.push_back(gridList[k]);
        k++;
    }
    return boundary_grid;
}

// template<unsigned int DIM>
// class UnstructedMesh: public Mesh<DIM>
// {
// public:
//     // virtual void readMesh();
//     // virtual void writeMesh();
//     // virtual void distributeDofs(int _order_of_element);
//     virtual void InputDofData() = 0;
//     virtual void InputEdgeData() = 0;
//     virtual void InputEleData() = 0;
// };

// class EasyMesh :public  UnstructedMesh<2>
// {
// private:
//     string _MeshName;
//     long int _ndofs;
//     long int _nele;
//     long int _nbndedge;
//     vector<Dofs<2>> _DofsList;
//     vector<Dofs<2>> _BndDofs;
//     vector<vector<int>> _ElementData;
//     vector<vector<int>> _BndEdges;

// public:
//     EasyMesh();
//     EasyMesh(string meshname);
//     void InputDofData();
//     void InputEdgeData();
//     void InputEleData();
//     int getNumGrid(){return _nele;};
//     int getTotalNumDofs(){return _ndofs;};
//     int getBndedgeNum(){return _nbndedge;};
//     Dofs<2> getDofs(int i){return _DofsList[i];};
//     vector<int> getElement(int i){return _ElementData[i];};
//     vector<Dofs<2>>& getBoundary(){return _BndDofs;};
//     //bool IsBoundary(int index){return !(std::find(_BndDofs.begin(),_BndDofs.end(),index) == _BndDofs.end());};
//     vector<int> getBndEdge(int i){return _BndEdges[i];};
// };

// EasyMesh::EasyMesh(string meshname)
// {
//     _MeshName = meshname;
//     InputDofData();
//     InputEdgeData();
//     InputEleData();
// }

// void EasyMesh::InputDofData()
// {
//     ifstream fin;
//     fin.open(_MeshName+".n");
//     string line;
//     stringstream ss;
//     getline(fin, line);
//     ss.clear();
//     ss.str(line);
//     ss >> _ndofs;
//     _DofsList.resize(_ndofs);

//     int a, bndflag;
//     double x, y;
//     for (int i = 0; i < _ndofs; i++)
//     {
// 	getline(fin, line);
// 	ss.clear();
// 	ss.str(line);
// 	ss >> a >> x >> y >> bndflag;

// 	Point<2> temp_point({x, y});
// 	Dofs<2> temp_dof(temp_point, i);
// 	_DofsList[i] = temp_dof;
// 	if (bndflag)
// 	    _BndDofs.push_back(temp_dof);
//     }
// }

// void EasyMesh::InputEdgeData()
// {
//     ifstream fin;
//     fin.open(_MeshName+".s");
//     string line;
//     stringstream ss;
//     getline(fin, line);
//     ss.clear();
//     ss.str(line);
//     int _nedge;
//     ss >> _nedge;
//     int a, d0, d1, bndflag;
//     for (int i = 0; i < _nedge; i++)
//     {
// 	getline(fin, line);
// 	ss.clear();
// 	ss.str(line);
// 	ss >> a >> d0 >> d1 >> a >> a >> bndflag;
// 	if (bndflag)
// 	{
// 	    vector<int> temp_edge(2);
// 	    temp_edge[0] = d0;
// 	    temp_edge[1] = d1;
// 	    _BndEdges.push_back(temp_edge);
//  	}
//     }
//     _nbndedge = _BndEdges.size();
// }

// void EasyMesh::InputEleData()
// {
//     ifstream fin;
//     fin.open(_MeshName+".e");
//     string line;
//     stringstream ss;
//     getline(fin, line);
//     ss.clear();
//     ss.str(line);
//     ss >> _nele;
//     _ElementData.resize(_nele);

//     int a, d0, d1 ,d2;
//     for (int i = 0; i < _nele; i++)
//     {
// 	getline(fin, line);
// 	ss.clear();
// 	ss.str(line);
// 	ss >> a >> d0 >> d1 >> d2;
// 	_ElementData[i].resize(3);
// 	_ElementData[i][0] = d0;
// 	_ElementData[i][1] = d1;
// 	_ElementData[i][2] = d2;
//     }
// }

#else
#endif
